% ESTIMA_LOOP:
% Function can also be called with one giant function handle, may be
% slower.
% estfunc = @(x)lmj_postmin(x, parstrval, parposest, ...
%                  partrmat_est, funcmod, Y, trainvec, prpar, prss, ...
%                  solveopt, addsol, ssposest, flag_transform);
%
% NUMPAR            Number of total parameters (including calibrated)

% 1. Possibly start a DIARYFILE
% =========================================================================
if flag_diary==1
    diaryfname=['log ',strdate(1)];
    % Create a Diary File of the Optimization
    delete(diaryfname);
    diary(diaryfname);
    disp('__________________________________________________');
    datestr(now);
end

% 2.Organize Outputs
% =========================================================================
logprior_mat(:,2)= -1e15*ones(nstrvals,1);
loglikel_mat(:,2)= -1e15*ones(nstrvals,1);
logpostd_mat(:,2)= -1e15*ones(nstrvals,1);
% PARMODEMAT     Matrix of output vectors
% OPTIMAT        Matrix with info on optimization
%                First columns number of iterations
%                Second column return flags
%                Third time in minutes
parmodemat = zeros(numpar, nstrvals);
optimat    = zeros(3, nstrvals);

% 2.
% Begin estimation Loop
% Will work with TRY and CATCH later on but not yet for debugging

% 3. Start a cell that can save output as we go along
% TAB_LOOP
% =========================================================================
tab_loop=emptycell(numpar+7,nstrvals+1);
tab_loop(1:numpar,1)=parnames;
tab_loop(numpar+2:end,1)={'logPost';'logLikel';'Number';'Iteration count';'exitflag';'minutes'};
cd(outpath);save tab_loop;cd(cucd);


% 4. Assign default values to ESTIMOPT if they do not exist
% =========================================================================
if ~exist('estimopt', 'var'); estimopt = []; end;
[junk,estimopt]=ch_field(estimopt,'TolFun',1e-5);
[junk,estimopt]=ch_field(estimopt,'MaxFuneval',200000);
[junk,estimopt]=ch_field(estimopt,'MaxIter',400);
[junk,estimopt]=ch_field(estimopt,'Display','on');
if case_estima == 4; % when using Simann
    % pass in special inputs
    [junk,estimopt]=ch_field(estimopt,'sa_t',1e7);
    [junk,estimopt]=ch_field(estimopt,'sa_rt',0.9);
    [junk,estimopt]=ch_field(estimopt,'sa_nt',5);
    [junk,estimopt]=ch_field(estimopt,'sa_ns',20);
    [junk,estimopt]=ch_field(estimopt,'rseed',999);
    [junk,estimopt]=ch_field(estimopt,'sa_neps', 4);
    [junk,estimopt]=ch_field(estimopt,'sa_eps', 1e-5); % convergenece criteria
    [junk,estimopt]=ch_field(estimopt,'sa_maxeval', 12*1e6);
    [junk,estimopt]=ch_field(estimopt,'sa_maxiter',  3*1e6); %3*1e6);
end

% 5. Assign default values to a PATH in case connection breaks down
% ------------------------------------------------------------------------
if isunix==0
    backup_path=cr_dir('C:','matlabrepo',root,spec,subf);
else
    backup_path=outpath;
end
% =========================================================================
% 6.
% Preliminary assigments and truncations before loop begins
% Formerly known as INDMAT, this is the matrix of transformations
% PARVECMODE is an empty vector with the calibrated values in
% =========================================================================
partrmat_est          =partrmat(parposest, :);
parvecmode            =zeros(numpar, 1);
parvecmode(parposcal) =parcalval;
errcount              =0;

tiempo=0; 

for ii=1:nstrvals;
    
    disp('===================================================');
    disp('===================================================');
    disp(['Starting Candidate Value number : ',num2str(ii) ] );
    xzero=parstrvals(parposest,ii);
    flag_optimok=1;
    
    tic;
    estima_loop;
    time=toc/60;
    tiempo=tiempo+time; 
    
    if flag_optimok==0
        continue
    end
    
    % =========================================================================
    % 8. OUTPUT
    % Output that must be produced by each code
    % -----------------------------------------
    % XESTMODE:  Estimated value
    % ESTIMITER: Number of Iterations
    % EXITFLAG : How routine ended
    % LPOSTMIN : Exit value of posterior (min)
    % =========================================================================
    
    parvecmode(parposest) = xestmode;
    parestmode = xestmode;
    
    cd(cucd); 
    
    % =================================================================
    % 8.1. Evaluate Posterior and Likelihood at the MODE
    %      Write them to the temporary TAB matrix
    % =================================================================
    [lpostd, likel] = lmj_postmax(xestmode,parvecmode,parposest,funcmod,Y,trainvec, prpar, prss, solveopt, addsol, ssposest);
    
    logpostd_mat(ii, 2) = lpostd;
    loglikel_mat(ii, 2) = likel;
    logprior_mat(ii, 2) = lpostd - likel;
    
    % PARMODEMAT     Matrix of output vectors
    parmodemat(:, ii) = parvecmode;
    
    optimat(1, ii) = estimiter;
    optimat(2, ii) = exitflag;
    optimat(3, ii) = time;
    
    % =================================================================
    % 8.2.
    % Print to screen
    % ===================================================================
    try
        disp(['Starting parameter ', num2str(ii), ': Estimated values are '])
        printcell([parnames(parposest), num2cprec(parestmode)])
        printcell([{'logPost'; 'logLht'}, num2cprec([lpostd; likel])])
    catch
        disp('errors in DISPLYAING RESULT')
    end
    
    % =================================================================
    % 8.3.
    % Write TAB_LOOP
    % ===================================================================
    tab_loop(1:numpar,ii+1)=num2cprec(parmodemat(:,ii),15);
    tab_loop{numpar+2,ii+1}=num2str(logpostd_mat(ii,2));
    tab_loop{numpar+3,ii+1}=num2str(loglikel_mat(ii,2));
    tab_loop{numpar+4,ii+1}=num2str(ii);
    tab_loop{numpar+5,ii+1}= num2str(estimiter); % Iteration count
    tab_loop{numpar+6, ii+1}= num2str(exitflag); % Exit flag
    tab_loop{numpar+7, ii+1}= num2str(time);
    
    % =================================================================
    % 8.4.
    % Attemp to save output
    % ===================================================================
    try
        cd(outpath);
        save mat_loop logprior_mat logpostd_mat logprior_mat ...
            parmodemat optimat loglikel_mat
        if ~isunix
            xlswrite('tab_loop', tab_loop);
        end
        cd(cucd);
    catch
        cd(backup_path)
        if ~isunix
            xlswrite('tab_loop', tab_loop);
        end
        save mat_loop logprior_mat logpostd_mat logprior_mat ...
            parmodemat optimat loglikel_mat
        cd(cucd)
        disp('Could not save the data')
    end
    
    
    % =====================================================================
end

disp(' ');
disp('____________________________________');
disp('End all draws');

dispaj('Total time in minutes=',sum(optimat(3,:)));

if flag_diary==1
    diary off;
    copyfile(diaryfname,outpath);
    cd(cucd);
    eval(['delete ',diaryfname,'.txt']); 
end
dispaj('Total Running time=',tiempo); 
cd(outpath); 
save workspace; 
cd(cucd); 

